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ABSTRACT 

We have determined the cosmological evolution of the density of active galactic nuclei (AGN) and of 
their Njj distribution as a function of the un-absorbed 2-10 keV luminosity up to redshift 4. We used 
the HELLAS2XMM sample combined with other published catalogs, yielding a total of 508 AGN. 
Our best fit is obtained with a luminosity-dependent density evolution (LDDE) model where low 
luminosity (L x ~10 43 erg s" 1 ) AGN peak at z~0.7, while high luminosity AGN (L x >10 45 erg s" 1 ) 
peak at z~2.0. A pure luminosity evolution model (PLE) can instead be rejected. 

There is evidence that the fraction of absorbed (Nn>10 22 cm" 2 ) AGN decreases with the intrinsic 
X-ray luminosity, and increases with the redshift. 

Our best fit solution provides a good fit to the observed counts, the cosmic X-ray background, and 
to the observed fraction of absorbed AGN as a function of the flux in the 10 _15 <S , 2-io<10 -10 erg 
cm" 2 strange. We find that the absorbed, high luminosity (Lx> 10 44 erg s _1 ) AGN have a density 
of 267 deg -2 at fluxes S2-io>10 -15 erg cm" 2 s _1 . Using these results, we estimate a density of 
supermassive black holes in the local Universe of pbh = 3.2 h 2 x 10 5 M Q Mpc -3 , which is consistent 
with the recent measurements of the black hole mass function in the local galaxies. 

Subject headings: diffuse radiation — galaxies: active — galaxies: evolution — quasars: general — 
surveys — X-rays: diffuse background 



1. INTRODUCTION 

The understanding of the history of accretion in the 
Universe and of the formation of massive black holes and 
their host galaxies relies on the measurement of the active 
galactic nuclei (AGN) space density and evolution. 

According to the AGN unified model (Antonucci 1993) 
the viewing angle between the observer and the symme- 
try axis of the nuclear structure is responsible for the 
different classification. In type 1 AGN the central en- 
gine is directly visible. Both the broad and narrow line 
emitting regions are detected in the optical spectra along 
with a soft un-absorbed X-ray spectrum. On the con- 
trary, a type 2 AGN classification arises when the broad 
line region and the soft X-rays are obscured by a dusty 
torus. 

Until a few years ago the best measurements of the 
cosmological evolution of the AGN luminosity function 
were essentially limited to optically (e.g. La Franca & 
Cristiani 1997, Croom et al. 2004), and soft X-rays (e.g. 
Maccacaro et al. 1991, Miyaji et al. 2000) selected type 
1 AGN. While there is evidence that type 2 AGN are 
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about a factor four more numerous than type 1 AGN (e.g. 
Maiolino & Rieke 1995; Risaliti et al. 1999), their rela- 
tive space density beyond the local Universe is basically 
unknown. Assuming that the cosmological evolution of 
type 1 and 2 AGN is the same, it was possible to simulta- 
neously reproduce the X-ray background spectrum and 
the X-ray counts (e.g. Setti & Woltjer 1989; Comastri et 
al. 1995). This simple picture was later slightly modified 
in models where the fraction of type 2 AGN was assumed 
to increase towards higher redshifts (e.g. Pompilio et al. 
2000, Gilli et al. 1999, 2001). The selection of complete 
samples of type 2 AGN is a difficult task. In the opti- 
cal they are often so dim that only the light of the host 
galaxy is visible; at z >1 even the latter has usually R > 
24. In the soft X-rays bands even hydrogen column den- 
sities, Nn, of the order of io 21 " 22 cm" 2 may strongly 
suppress the flux. In the hard (2-10 keV) X-rays type 2 
AGN selection is less biased against, though the absorp- 
tion due to large Ah column densities (10 23-24 cm" 2 ) is 
not negligible especially at low redshifts. 

Early attempts to compute the hard X-ray luminosity 
function, based on ASCA and Beppo-SAX observations 
(Boyle et al., 1998; La Franca et al. 2002 respectively) 
indicated a strong evolution for type 1 AGN, with a rate 
similar to that measured in the soft X-rays. Unfortu- 
nately the low spatial resolution of the X-ray detectors 
prevented an unambiguous identification of the type 2 
AGN optical counterparts, thus hampering a reliable de- 
termination of the type 2 AGN space density. 

Thanks to the high sensitivity and spatial resolution 
of the hard X-ray detectors on board XMM- Newton and 
Chandra, it has become possible to carry out AGN sur- 
veys less biased against X-ray absorption and with more 
secure optical identifications. 
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However, already at fluxes fainter than S2-io~10 
erg cm -2 s _1 , a sizeable fraction of the X-ray sources 
have optical magnitudes fainter than the spectroscopic 
limit of 8-10 meter class optical telescopes, and thus the 
measure of their distance has to rely on photometric red- 
shifts, when it is not impossible altogether. 

For these reasons, although the Chandra Deep Field 
North (CDF-N; Alexander et al. 2003) and the Chan- 
dra Deep Field South (CDF-S; Giacconi et al. 2002) sur- 
veys have resolved a fraction of the 2-10 keV XRB as 
large as 85-90% (see also Brandt and Hasinger 2005), 
a clear picture of the AGN evolution able to reproduce 
the whole set of observational constraints (i.e. soft and 
hard X-ray counts, X-ray background, Ah and redshift 
distributions) is still missing. 

Attempts to take into account the redshift incomplete- 
ness of X-ray selected AGN have been carried out by 
Cowie et al. (2003), Fiore et al. (2003), Barger et al. 
(2005) combining data from deep and shallow surveys. 
They independently demonstrated that the AGN number 
density for luminosities lower than ~10 44 erg s _1 peaks 
at a lower redshift than that of high luminosity objects. 
Making use of an almost complete sample of 247 AGN 
from Chandra, ASCA and HEAOl surveys above a lim- 
iting flux of S 2 -io>3.8 x 10~ 15 erg cm~ 2 s _1 Ueda et 
al. (2003) were able to estimate the hard X-ray luminos- 
ity function (HXLF) up to z = 3. They found that the 
fraction of the X-ray absorbed AGNs decreases with the 
intrinsic luminosity and that the evolution of the AGN 
HXLF is best described by a luminosity-dependent den- 
sity evolution (LDDE). Very similar results were also ob- 
tained by Hasinger et al. (2005) using an almost complete 
sample of soft X-ray selected type 1 AGN. 

In this paper we expand the study carried out by Fiore 
et al. (2003) with the aim to compute the shape and 
evolution of the HXLF and Ah distribution of all AGN 
with Ah<10 25 cm -2 up to z ~ 4. To reach such a goal it 
is necessary to cover the widest possible range in the Lx- 
z-Ah space, and to take into account all possible selection 
effects. For these reasons we have used a large AGN 
sample (about 500 objects) four times deeper than the 
Ueda et al. (2003) sample. A new method to correct for 
the spectroscopic incompleteness of faint X-ray sources 
is presented and discussed in detail. The selection effects 
due to X-ray absorption are also specifically discussed 
and estimated by an appropriate X-ray "K-correction" 
term. 

The paper is structured as follows: in Section 2 we 
describe the adopted X-ray samples; in Section 3 the 
method to compute the HXLF is discussed. The results 
are presented in Section 4, discussed in Section 5 and 
summarized in the last section. 

Throughout this paper we call AGN all objects with 
an intrinsic (corrected for Ah absorption) 2-10 keV X- 
ray luminosity larger than 10 erg s _1 . In the last few 
years evidence for a mismatch between optical (Type 
1/2) and X-ray (un-absorbed/absorbed) classification 
has emerged (e.g. Fiore et al. 2000). In this paper 
we refer to AGN1 if broad emission lines (rest frame 
FWHM>2000 km s" 1 ) are present, while all remaining 
objects (with or without narrow emission lines in the 
optical spectrum) are called AGN2. If the rest frame 
column density is larger than 10 22 cm -2 the AGN is 
classified as absorbed. The adopted limit is well above 



TABLE 1 

The samples 



Sample 


Flux limit 


No 




1 ' I, i ni 


(1) 


(2) 


(3) 


(4) 


(5) 


HEAOl 


2.9X10" 11 


31 


31 




AMSSn 


3.0X1CT 13 


71 


74 




HBS28 


2.2xl0" 13 


27 


27 




H2XMM a 


8.0X1CT 15 


120 


103 ( 93) 


23.65 


H2XMM b 


8.0X1CT 15 


110 


44 


21.40 


Lockman 


2.6xl0" 15 


55 


41 ( 39) 


23.50 


CDFN 


l.OxlO" 15 


146 


108 (102) 


24.65 


CDFS 


l.OxlO" 15 


127 


102 ( 98) 


25.00 


Note. — 


In column (2) 


wc give 


the fiux limit 


of the 



samples in units erg cm" 2 s — 1 . In column (3) ve give 
the total number of sources. In column (4) we give 
the number of sources brighter than the spectroscopic 
limit (in parenthesis those having redshift). In column 
(5) we give the spectroscopic completeness magnitude. 

a ldf sample (Fiore et al. 2003). 
b 0.5df sample (Cocchia et al. 2005). 



the typical X-ray absorption by host galaxy gas (disk, 
starburst regions, etc.) thus ensuring that the measured 
column is most likely related to nuclear obscuration. Un- 
less otherwise stated, all quoted errors are at the 68% 
confidence level. We assume Hq — 70 km s _1 Mpc -1 , 
O m = 0.3 and Q\ = 0.7. 



2. SAMPLES 

In order to cover the widest possible range of lumi- 
nosities and redshifts we combined the HELLAS2XMM 
sample (Fiore et al. 2003) with other existing flux limited 
samples which allowed the estimates of the rest frame Ah 
column density of each AGN. Whenever possible, the col- 
umn density and the photon index (r) were determined 
with a proper spectral analysis. Otherwise, we assumed 
r = 1.8, and used the hardness ratio to measure the z=0 
column density (Nho, see also the discussion about the 
uncertainties of this approach in §4.1.1). The rest frame 
column density (Ah) was then estimated by the relation 
Log(N H ) = Log(N H o) + 2A2Log(l + z), which makes 
use of the Morrison & McCammon (1983) cross sections, 
including also the effects of the absorption edges, and as- 
sumes solar abundances from Anders & Grevesse (1989). 

For those samples whose optical spectroscopic iden- 
tifications are incomplete, we chose the optical magni- 
tude limit at which the samples are almost spectroscop- 
ically complete. The incompleteness is 6% in the HEL- 
LAS2XMM, Lockman, CDF-N and CDF-S samples). In 
these cases (as the X-ray-optical flux distribution of the 
sources without redshift is almost similar to that of the 
spectroscopically identified sources, and the fraction of 
the unidentified sources is small) the sky coverage has 
been reduced according to the fraction of spectra avail- 
able. Table ^ contains a summary of the characteristics 
of each sample. The distribution in the Lx~z space of all 
AGN from the spectroscopically complete sub-samples 
used in our analysis are shown in Figure ^ while Figure 
shows their distribution in the S^-R plane. 
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log L x (2-10 keV) (deabsorbed) 

Fig. 1. — Lyi—z plane for all AGN used in this analysis. Different 
symbols corresponds to different surveys, as labeled in the top left 
corner. Absorbed sources are also highlighted by a cross. 
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2.1. T/ie HELLAS2XMM sample 

We used the HELLAS2XMM ldf (1 degree field) sam- 
ple (Fiore et al. 2003) plus the recently available ex- 
tension of 0.5 deg 2 (HELLAS2XMM 0.5df Cocchia ct 
al. 2005). The HELLAS2XMM ldf sample contains 122 
sources, serendipitously detected in five XMM-Newton 
fields with S* x (2-10 keV)> 0.8 x 10~ 14 erg cur 2 s" 1 . In 
our analysis we used the fluxes and the column densi- 
ties measured by X-ray spectral analysis (Perola et al. 
2004). Among the 122 sources we discarded one star 
(object n. 0537006) and one extended source (object n. 
26900013). For three sources with low signal-to- noise the 
hardness ratio and redshift were used to estimate the rest 
frame Ah- In summary, the sample contains 120 sources, 
115 optically identified, and 95 with measured redshift 
and optically classified. We restricted our analysis to the 
sources brighter than R=23.65. Down to this limit 93 
out of 103 sources have been spectroscopically identified. 

The HELLAS2XMM 0.5df sample consists of 110 ob- 
jects brighter than Sx(2-10 keV)= 8 x 10~ 15 erg cm~ 2 
s _1 . Among them, 44 sources brighter than R=21.4 (but 
otherwise randomly selected) have been spectroscopically 
identified. 

2.2. The Piccinotti sample 

The Piccinotti sample is the brightest included in our 
analysis. It has been obtained through observations car- 
ried out by the HEAOl satellite, and contains 31 sources 
selected over an area of 26919 deg 2 down to Sx(2-10 
keV)= 2.9 x 10~ n erg cm" 2 s" 1 (Piccinotti et al. 1982). 
The column densities have been taken from the litera- 
ture, and are derived from X-ray spectral analyses. 



Fig. 2. — R-band magnitude versus the 2-10 keV X— ray flux 
for all sources of the samples used in this analysis. The dashed 
lines are the spectroscopic limits of completeness adopted in our 
analysis. 



2.3. The AMSSn sample 

The AMSSn sample consists of 74 AGN at fluxes 
brighter than S x {2-10 keV)= 3 x 10~ 13 erg cm~ 2 
s _1 (Akiyama et al. 2003). The total area covered is 
45 deg 2 at the fainter fluxes and rises up to ^69 deg 2 
at bright fluxes. The Ah column densities have been 
derived from the hardness ratios values. 

2.4. The HBS28 sample 

The HBS28 sample (Caccianiga et al. 2004) consists 
of 27 AGN and 1 star selected in the 4.5-7.5 keV band. 
The sources are brighter than S x (2-10 keV)= 2.2 x 10~ 13 
erg cm~ 2 s _1 (assuming r=1.8) and have been selected 
over 82 XMM-iVewton pointed fields, corresponding to a 
total of 9.756 deg 2 . All sources have been spectroscop- 
ically identified, and their column densities have been 
measured through X-ray spectral fits. 

2.5. The Lockman Hole sample 

The Lockman Hole sample consists of 55 sources se- 
lected within the 12 arcmin radius of the XMM-Newton 
observation. The sources are brighter than S x (2-10 
keV)= 2.6 x 10" 15 erg cm" 2 s'^Baldi et al. 2002). 
Optical identifications and X-ray spectral fits are from 
Mainicri et al. (2002). Spectroscopic redshifts and clas- 
sifications have been obtained for 41 objects, while 3 
sources have photometric redshifts. We restricted our 
analysis to the sources brighter than R=23.50. Down to 
this limit 39 out of 41 sources have been spectroscopically 
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identified. 



2.6. The CDF-N sample 

In order to reach almost spectroscopic completeness we 
have selected an X-ray bright subsample in the CDF-N. 
The subsample consists of 146 sources (see Table 1) se- 
lected within the 10 arcmin radius of the Chandra ob- 
servation (Alexander et al. 2003). The sky coverage 
reaches Sx(2-10 keV)> 10~ 15 erg cm~ 2 s^ 1 in the in- 
ner 5.85 arcmin radius, Sx(2-10 keV)> 2.49 x 10" 15 
erg cm~ 2 s _1 in the annulus between 5.85 and 6.5 ar- 
cmin radii, and S , x (2-10 keV)> 3.61 x 10~ 15 erg cm" 2 
s _1 in the annulus between 6.5 and 10.0 arcmin radii. 
We used both spectroscopic and spectro-photometric 
identifications and redshifts available from the literature 
(Barger et al. 2003). We restricted our analysis to sources 
brighter than R=24.65. Down to this limit 102 out of 108 
sources have been spectroscopically identified. The Ah 
column densities have been derived from the hardness 
ratios. 



2.7. The CDF-S sample 

Altough the CDF-S has been observed for 1 Ms in- 
stead of the 2 Ms spent in the CDF-N, we selected a 
spectroscopically complete X-ray bright subsample with 
the same sky coverage as for the CDF-N. Indeed, at our 
adopted flux limits, the difference in the exposure time 
does not affect the sky coverage. The sample consists 
of 127 sources (see Table 1; Giacconi et al. 2002 and 
Alexander et al. 2003). We used both spectroscopic and 
spectro-photometric redshifts available from the litera- 
ture (Szokoly et al. 2004; Zheng et al. 2004). Moreover, 
given that both Szokoly et al. (2004) and Zheng et al. 
(2004) identifications are based on the X-ray source cat- 
alogue of Giacconi et al. (2002), we have revised some 
optical/X-ray associations according to the improved as- 
trometry provided by Alexander et al. (2003). We re- 
stricted our analysis to sources brighter than R=25.00. 
Down to this limit 98 out of 102 sources have been spec- 
troscopically identified. The Ah column densities have 
been derived from the hardness ratios. 



3. METHOD 

We searched for a functional fit to the density of the 
AGN as a function of the un-absorbed 2-10 keV lumi- 
nosity (Lx), the rest frame absorbing column density 
(Ah), and the redshift (z). The method is based on the 
comparison, through \ 2 estimators, of the observed and 
expected numbers of AGN (in the Lx-z space) and of 
the Ah distributions, obtained from computations which 
take into account all the observational selection effects of 
the samples. 

Once a HXLF evolution model is assumed, the number 
of expected AGNs (E) in a given bin of the Lx-z-N^i 
space is the result of the sum, over the number of samples 
A samp , of the expected number of AGN in each sample 
taking into account the area coverage of each ith sample 
rti(L, Nh, z), the Ah distribution f(Lx,z;Nn), and a 
completeness function g(Lx, z, Nh, Ri), where Ri is the 
spectroscopic limit of completeness of the ith sample: 



$(Lx,z)f(Lx,z;N u ) x 



dV 



g(L x , z, N H , Ri)^i{L, N H , z)—dLogLxdzdN H . (1) 

3.1. The shape of the Luminosity Function 

In order to describe the evolution of the AGN, we used 
standard functional forms, such as the pure luminosity 
evolution (PLE) model and a luminosity-dependent den- 
sity evolution (LDDE) model (see next Section and, e.g., 
Boyle et al. 1998; Miyaji et al. 2000; La Franca et al. 
2002; Ueda et al. 2003). The HXLF, representing the 
number density per unit comoving volume and per unit 
Log Lx, as a function of Lx and z, was expressed as: 



d<P{Lx,z) 
dLogLx 



(2) 



We adopted a smoothly-connected two power-law form 
to describe the present-day HXLF, 



d$(L x ,z = 0) 
dLogLx 



= AKLx/L^ 1 + (Lx/L^ 2 }- 1 - (3) 



3.2. The K-correction 

In order to convert the observed 2-10 keV fluxes (Sx) 
to the intrinsic 2-10 keV luminosities {Lx) and vice- 
versa, for each observed or "expected" AGN with a given 
Ah, a K-correction has been computed by assuming a 
photon index T = 1.8, an exponential cutoff (e~ E / Ec ) at 
Ec = 200 keV, and the corresponding photoelectric ab- 
sorption (see §4.1.1 for a discussion on the use of different 
K-corrections) . 

3.3. The completeness function 

All the faint samples used in our analysis (HEL- 
LAS2XMM, Lockmann, CDF-S, CDF-N) are nearly 
spectroscopically complete down to a certain optical limit 
magnitude (72=21.4 — 25, see TableQ]). In order to com- 
pute the number of expected AGN in a certain bin of the 
Lx-z-Nn space, we introduced the completeness func- 
tion g(Lx, z, Nh, R) which provides the probability that 
a given AGN with luminosity Lx, redshift z and column 
density Ah , had an apparent R-band magnitude brighter 
than the spectroscopic limits of completeness R of each 
sample. 

For this reason we derived an empirical relationship 
between the un-absorbed X-ray luminosity Lx and the 
optical luminosity Lr 9 for AGN1 and AGN2, and mea- 
sured their spread (see Figure For AGN1 we found: 

LogL R = 0.959(±.025) x LogL x + 2.2(±1.1), (4) 

with a la dispersion of 0.48 (in LogL# units) around 
the best fit solution. The linear correlation coefficient is 

9 The Lr luminosity is in erg s — 1 (vL v ), computed at 660 nm, 
where the flux is /[erg s" 1 cm" 2 Hz" 1 ] = 2.84 X 10~ 20 X l0-° iR 
(Zombeck 1990). 
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Fig. 3. — LogLx-LogL H relation for optical AGN1 and AGN2. 
The continuous lines correspond to eq. l41andl^l 



r=0.773, corresponding to a negligible (<10~ 13 ) proba- 
bility that the data are consistent with the null hypoth- 
esis of zero correlation. For AGN2 a flatter relation was 
found: 



LogL R = 0.462(±.026) x LogL x + 23.7(±1.1) 



(5) 



with a ler dispersion of 0.40 (in LogL^ units), and a lin- 
ear correlation coefficient r=0.462, again corresponding 
to a negligible (< 2xl0~ 13 ) probability that the data are 
consistent with the null hypothesis of zero correlation. In 
order to compute the above relationships a linear least 
squared method with errors (assumed 0.2 dex) in both 
axes has been used. The difference between the two rela- 
tions should be attributed to the dominance in the opti- 
cal of the AGN component in the AGN1, which produces 
an almost linear relationship between X-ray and optical 
luminosity (see La Franca et a!., 1995 for similar results 
in the soft X-rays). In AGN2, where the nucleus is ob- 
scured, the optical luminosity is instead dominated by 
the host galaxy (see also Fiore et al. 2003). 

For each pair of un-absorbed X-ray luminosity and red- 
shift, the above relationships (with their spreads) can be 
used to compute the probability of an AGN to appear 
brighter than a certain optical magnitude, and thus be 
spectroscopically identified. The observed spreads of the 
two relationships are due to a combination of the intrin- 
sic spread with the observational uncertainties. Given 
our aims, both effects should be taken into account, and 
we have thus not subtracted the contribution of the ob- 
servational uncertainties from the spread estimates. To 
choose which Lx-Lr relationship to use (eq. 2] or ec l- 
EJ, we need also to know the probability of an AGN to 
appear as an AGNI (or, its complement, an AGN2) as a 
function of Lx, Ah and z: QI(Lx, z, Ah). This proba- 
bility was estimated from the sample itself as described 
below. 

Figure 0| shows the distribution of the observer frame 
column density N#o as a function of Lx for AGNI and 




42 43 44 45 100 50 

Log Lx (2-10 keV) [erg s" 1 ] CDF 

Fig. 4. — LogLx-LogAf^fo plane at z < 1 (bottom), 1 < z < 
2 (middle), z > 2 (up). Open squares are optical AGNI, filled 
circles are optical AGN2. On the right side the dashed lines are the 
cumulative distribution functions of the Nj^q values for AGNI with 
Lx>10 43 ' 5 erg s —1 , while the continuous lines are the complement 
of the cumulative distribution function of the N^-q values for AGN2 
with L x >10 43 - 5 erg s^ 1 . 



AGN2, in three redshift intervals. Here we do not use the 
rest frame Ah, but instead the observer frame N# which 
is equivalent to an hardness ratio (see also Hasinger 
2003). As can be seen in Figure 03 the probability to 
find an AGNI is not only dependent on Nho, but de- 
pends also on the luminosity. The probability to find an 
AGNI increases with increasing luminosities, and there 
is a relevant fraction of low luminosity (Lx<10 43 erg s" 1 ) 
un-absorbed objects which are AGN2, while a fraction of 
the high luminosity (L x >10 45 erg s _1 ) absorbed objects 
are AGNI. This result, if it is not due to the contamina- 
tion by the galaxy light in the lower luminosity AGN2, is 
against the simplest version of the AGN unified model. 
The analysis of this issue is beyond the scope of this pa- 
per (see Panessa & Bassani 2002, Page et al. 2003, Steffen 
et al. 2003, Ueda et al. 2003, Brusa et al. 2003, Perola 
et al. 2004 and Barger et al. 2005 for similar results and 
discussions. See also §4.6) 

As Figure0]shows, there is no evidence of a dependence 
on redshift of the distribution of AGNI and AGN2 as a 
function of Lx and Nho- We have thus estimated the 
probability of an AGN to appear as an AGNI as a func- 
tion of Lx and Nho only, by assuming no dependence on 
redshift. This probability has been estimated as a func- 
tion of L x in two bins of N ff0 10 : at N H0 <10 21 - 5 cm -2 

10 We chose to use, here, the observed column densities (Nho) 
instead of the intrinsic ones (N-r) in order to eliminate the de- 
pendencies on the redshift. A constant (with z) N//o=10 21 ' 5 
cm -2 separation limit corresponds to a shifts of this limit to- 
wards higher values of JVh with increasing redshift (as the in- 
trinsic and the z=0 column densities are related by the equation 
log(N H ) = log(N H0 ) + 2A2log(l + z)). We will come back to this 
point in the next Sections. However, we wish to stress here that the 
above relationships have been derived only in order to correct the 
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43 44 
Log L x (2-10 keV) [erg s" 1 ] 

Fig. 5.— Fraction of optical type 1 AGN with N H0 >10 21 ' 5 cm -2 
(up), and N£fo<10 21,5 cm -2 (bottom) as function of the intrinsic 
luminosity Lx- 



in 0.4 




43 44 45 46 0.5 1 
Log Lx (2-10 keV) z 

Fig. 6. — Observed fraction of absorbed (iV H >10 22 cm -2 ) AGN 
as a function of Lx (left panel) and z (right panel). The dotted 
lines are the intrinsic fractions assuming a constant flat Nn distri- 
bution in the range 10 20 <-/Vjj<10 25 cm -2 . The dashed lines are 
the expectations taking into account the selection effects. 



and Njy >10 21 ' 5 cm -2 . The values of the probability of 
an AGN to appear as an AGN1 in these two Nho inter- 
vals are shown in Figure El 

We caution the reader that, due to inhomogeneities on 
the quality of the spectroscopic classification of the sam- 
ples used, the above measure of the fraction of AGN1 as 
a function of Lx and N#o has uncertainties that are diffi- 
cult to quantify. However, these estimates are only used 
to derive which fractions of the unidentified AGN will 
follow the two Lx-Lr relationships shown in eq. 0] and 
El The absence of many outliers in the Lx- Lr relation- 
ships for AGN1 and AGN2 shown in Figure El demon- 
strates qualitatively that classification errors should not 
be very large. This, in principle, does not imply that 
more accurate spectroscopy would not change the opti- 
cal classification of the AGN, but that the spectroscopy 
is accurate enough, for our purposes, to decide which 
of the two Lx- L-r relationships the AGN would follow. 
However, the completeness correction is computed under 
the assumption that the measured fraction of AGN1 as a 
function of Lx and N#o and the derived two Lx- Lr re- 
lationships for AGN1 and AGN2 hold also for the higher 
redshift, optically fainter, unidentified population. We 
will discuss in §4 how much the uncertainties on these 
assumptions would affect our results. 

In summary, the completeness function 
g(Lx,z,NH,Ri) was computed as follows: for 
each given AGN having intrinsic luminosity Lx, 
redshift z, and absorption column-density Ah, a) 
the N#o was derived according to the equation 
Log{N m ) = Log(N H ) - 2A2Log{\ + z), b) the prob- 
abilities to be an AGN1 (Ql) and AGN2 (1-Q1) were 
estimated according to the values plotted in Figure El 
and then c) according to eq. El and El and their spreads, 
the two probabilities (for the fraction of AGN1 and 
AGN2) to be brighter than the spectroscopic limit R^ of 
the ith sample were computed and summed. 

samples for spectroscopic incompleteness, and that Nj^q = 10 21 ' 5 
cm -2 should not be meant as the working separation limit be- 
tween absorbed and un-absorbed AGN, which (as defined in §1) 
instead is N ir =10 22 cm~ 2 . 



3.4. The Ah function 

To describe the distribution of the spectral parame- 
ters of the AGNs at a given luminosity and redshift, we 
introduced the Ah function, f(Lx, z; Ah), a probability- 
distribution function for the absorption column-density 
as a function of Lx and z. The Ah function (in LogAn -1 
units) is normalized to unity at each redshift, over the 
A H interval 20<LogA H <25: 




/(Lx^;A H )dLogA H = 1. (6) 



The objects have been grouped into 5 bins of Ah, 
ALogAn=l wide, and centered at LogAn= 20.5, 21.5, 
22.5, 23.5, 24.5. The first bin includes all the AGNs hav- 
ing A H <10 21 cm" 2 . 

In Figure El the observed fraction of absorbed 
(Ah>10 22 cm -2 ) AGN as a function of Lx and z is 
shown. The dotted lines correspond to the fraction of 
absorbed objects if a flat Ah distribution in the range 
10 20 <N#<10 25 cm -2 were assumed, with no selection 
effects taken into account. The dashed lines show our 
predictions when these effects are included. Such a model 
does not provide a good fit to the data points, where a 
decrease with the intrinsic luminosity and an increase 
with the redshift is observed. This behavior is also ev- 
ident in Figure d Bearing in mind that the Ah esti- 
mates are affected by uncertainties that can be as large 
as one decade, from the analysis of Figure it appears 
that the assumption of a flat Ah distribution produces 
an expected distribution roughly in agreement with the 
observed one at Ah>10 21 cm -2 . Hence, the observed 
change of the fraction of absorbed AGN as a function of 
Lx and z (see Figure EJ, could be mainly attributed to a 
change of the fraction of AGN with Ah<10 21 cm~ 2 . This 
will be our working hypothesis, which we will analyze in 
the next Sections. 

We have thus assumed a flat Ah distribution between 
Ah= 10 21 cm~ 2 and Ah=10 25 cm' 2 , while allowing the 
fraction of objects with Ah<10 21 cm~ 2 to vary. We in- 
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logN H [cm-*] 



Fig. 7. — Ah distributions in various luminosity and rcdshift 
intervals. Top. High luminosity and low luminosity AGN at z<l. 
Bottom. High and low redshift AGN. The continuous lines are the 
observed distributions, the dotted lines are the assumed (constant 
flat) Ah distributions, while the dashed lines are the expectations 
taking into account the selection effects. 



traduced a linear dependence of the fraction of objects 
with A H <10 21 cm- 2 (* = f{L x ,z;LogN H < 21)) on 
both LogLx and z: 



^(L x ,z) = ^[{LogL x -U)/3 L + l}[(z-0.5)f3 z + l], (7) 

where ip is the fraction of objects with Ah<10 21 cm -2 
at L x = 10 44 erg s" 1 and z=0.5, and (3l and (i z are the 
slopes of the linear dependences on L x and z respec- 
tively. This choice is the simplest possible according to 
the quality of the data. The function holds for the ranges 
0.25 < z < 2.75 and 10 42 - 5 <L X <10 45 5 erg s" 1 . At red- 
shifts and luminosities outside these ranges the fraction 
was kept constant, equal to the values assumed at the 
limits of the ranges. Obviously ^ could take all values in 
the range [0,1]. This corresponds to an allowed fraction 
of absorbed objects (Ah>10 22 cm -2 ) in the range 0% - 
75%. Indeed, according to eq. Eland the assumption of 
a flat Nh distribution at Ah>10 21 cm -2 , the fraction of 
absorbed AGN turns out to be: 



N(22 < LogN H < 25) 
N(LogN H < 25) 



0.75(1 - *) 



(8) 



As clear from Figure[7] no object with Ah>10 25 cm" 2 
is either expected or observed. Thus, we limited our 
statistical analysis of the evolution of the AGN to the 



objects having Ah<10 25 cm~ 2 (see eq. EJ. However, 
when we will predict the number counts, X-ray back- 
ground and the accretion history (§4.3, §4.4 and §5.2), we 
will include in the Ah distribution a number of objects 
with 10 25 <A r H<10 26 cm -2 , equal to that in the interval 
10 24 <A H <10 25 . 

3.5. x 2 fitting 

In order to find the best fitting model we choose two 
X 2 estimators as figure of merit functions. 

The first estimator (x 2 F ) is related to the shape and 
evolution of the HXLF and is obtained by comparing the 
expected and observed numbers of AGN in 24 bins, cov- 
ering the whole sampled Hubble space (L x -z; see Figure 
[S] for an example of the binning) . Computations have 
been carried out in the 0<z<4.5 redshift range, and in 
the 10 42 <L X <10 47 erg s _1 luminosity range. A total of 
508 AGN were used. 190 had the Ah column densities 
directly measured from X-ray spectroscopic analysis. 

The second estimator (xn h ) ^ s re l & ted to the Ah func- 
tion, f(L x ,z;Nn), i.e. the shape of the Ah distribution 
and its dependence on L x and z. One contingency table 
was created dividing the objects with column densities 
higher or lower than Ah= 10 22 cm -2 into 5 further bins 
in the L x -z space. The Xn h estimator was computed by 
comparing the expected and observed number of AGN 
in the total 10 ( 2x5 ) bins. 

The reasons for using two different x 2 estimators are: 
a) the number of objects is too small to construct a single 
X 2 estimator using bins in the three-dimensional space 
L x -z-Nn ; b) the two x 2 estimators cannot be summed 
as the data used are not independent. The shape and 
the evolution of the HXLF is only marginally dependent 
on the shape and evolution of the Ah distribution (we 
checked that the best fit parameters of the HXLF vary 
within the ler uncertainties when the parameters of the 
Ah distribution are left to vary within a 3c range of their 
best-fit values). 

The final fit was obtained by iteratively searching for 
the lowest values of Xlf anc ^ Xn h ^ n turn, until the 
changes on the two x 2 estimators were smaller then 0.1 11 . 
For each model the probabilities for Xlf an d Xn H ' ac ~ 
cording to the corresponding degrees of freedom, were 
computed. Confidence regions of each parameter have 
been obtained by computing A% 2 at a number of val- 
ues around the best fit solution, while letting the other 
parameters free to float (see Lampton et al. 1976). The 
68% confidence regions quoted correspond to Ax 2 =L0. 
Moreover, in order to use an un-binned goodness of fit 
test of the HXLF models, we used also a bi-dimensional 
Kolmogorov-Smirnov test (2D-KS, Fasano & Frances- 
chini 1987) on the Hubble (L x -z) space. 

4. RESULTS 

4.1. The LDDE model 
By introducing the evolution factor 



e(z) 



{1 + z)p 1 (z<z c ) 
e(z c )[(l + z)/(l + z c )F 2 (z>z c ), 



(9) 



11 This is a small enough interval, as the variance on the x 2 
estimator is 2N,j, where N<; are the degrees of freedom. Variations 
of A^ 2 =0.1 correspond to confidence levels of less than 2% and 3% 



for 



Xlf 



and 



Xn h respectively. 
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TABLE 2 
2-10 keV AGN LF parameters 



Model A a pi ^2 ~ t a LogL a b LogL* b 71 ~ 2 $ Tl Wz P lf (x 2 ) c P lf (KS) c Pn h (x 2 ) c XRB 2 _io d 



LDDE 



1 1.48 


4.37 


-1.19 


2.39 


0.20 


45.74 


44.26 


0.94 


2.35 


0.26 


0.00 


-0.00 


20 


6 


0.07 


1.76 


2 1.50 


4.39 


-1.14 


2.41 


0.20 


45.74 


44.25 


0.97 


2.36 


0.29 


0.43 


-0.00 


20 


6 


0.5 


1.78 


3 1.39 


4.48 


-1.19 


2.39 


0.20 


45.74 


44.26 


0.94 


2.35 


0.30 


0.00 


-0.33 


19 


5 


0.09 


1.75 


4 1.21 


4.62 


-1.15 


2.49 


0.20 


45.74 


44.25 


1.01 


2.38 


0.44 


0.62 


-0.51 


20 


7 


83 


1.81 


5 c 1.29 


4.85 


-1.03 


2.45 


0.22 


45.73 


44.23 


1.09 


2.44 


0.36 


0.67 


-0.00 


33 


21 


20 


2.16 


Errors 5% 


+0.26 
-0.26 


+0.72 
-0.68 


+0.82 
-0.68 


+0.04 
-0.03 


+0.58 
-0.63 


+0.18 
-0.18 


+0.08 
-0.10 


+0.13 
-0.11 


+0.04 
-0.05 


+0.14 
-0.13 


+0.14 
-0.17 










LDDE 


6 6.18 


3.22 




1.08 






43.79 


0.95 


2.74 


0.46 


0.64 


-0.58 


9 


3 


63 


2.54 


Errors 5% 


+0.13 
-0.26 




+0.08 
-0.06 






+0.15 
-0.11 


+0.06 
-0.07 


+0.27 
-0.23 


+0.04 
-0.05 


+0.14 
-0.13 


+0.14 
-0.17 











a In units of 10~ e h% Mpc~ 3 . 



In units of erg s 
c Probability values in % units. 

In units 10 — 11 erg cm -2 s _1 dcg~ 2 . 
e Only Piccinotti, AMSSn and CDF-N samples used. 



the pure density evolution (PDE) model is expressed as 



d$(L x ,z) = d$(L x ,0) 
dLogix dLogLx 



e(z). 



(10) 



The z c parameter represents the redshift at which the 
evolution stops, pi is the parameter characterizing the 
rate of the evolution, while p2 is usually negative and 
characterizes the rate of the counter-evolution of the 
HXLF at z>z c . 

The LDDE model is obtained by introducing in the 
PDE model a luminosity dependence of z Cl assumed to 
be a power law: 



z c (L x ) 



z*AL x /L a ) 



(L x > L a ) 
(L x < L a ). 



(11) 



The above parameterization has been introduced by 
Ueda et al. (2003) in order to allow for a change with 
luminosity of the redshift at which the density of AGNs 
peaks (see also Miyaji et al. 2000 for a similar LDDE 
parameterization). This behavior is also apparent in our 
data (see, e.g. Figure HJ. 

In order to plot the HXLF we adopted the 'W obs /A mdl 
method" (La Franca & Cristiani 1997), where the best- 
fit model multiplied by the ratio between the number 
of observed sources and that of the model prediction in 
each Lx~z bin is plotted. Although model dependent 
(especially when large bins are used), this technique is 
the most free from possible biases, compared with other 
methods such as the conventional 1/V a method. The at- 
tached errors are estimated from Poissonian fluctuations 
(lcr) in the observed number of sources according to the 
Gehrels (1986) formulae. 

We simultaneously fitted the parameters of the HXLF 
and of the possible dependencies of the Ah distribution 
on Lx and z. As shown in Table [3 the LDDE model 
provides a good fit to the data regardless of the adopted 
ATjj distribution (see Figures [S] and 0) ■ 

According to these fits, the redshift of the density peak 
of AGN increases with the luminosity, from z^0.5 at 



L x ~10 42 erg s" 1 up to z~2.5 at L x ~10 4b erg s" 1 . Out of 
the four proposed Ah distributions only fit # 4, provides 
a good fit to the whole data in the Lx-£-Ah space. The 
first model (fit # 1) searched for a constant value of 
the fraction of objects with Ah<10 21 cm~ 2 (^(Lx, z) = 
ip = constant, /3 L =/3 Z =0). The x 2 probabilities of the 
dependence of the Ah distributions on Lx and z reject, 
at more than 99.93% confidence level, this model. As can 
be seen in Figures |H]and[3 the data requires a decrease of 
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Fig. 8. — Density of AGN in luminosity bins as a function of 
redshift. The continuous lines are the best fit values of the LDDE 
model with an evolving Nn distribution depending on Lx and z 
(fit #4 in Table^J. Data have been plotted using the "AT° bs /jV mdl 
method" (see §4.1). The dashed areas are the largest allowed re- 
gions due to uncertainties in the completeness correction method 
used (see §4.1.1). 
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Fig. 9. — Density of AGN as a function of luminosity in four red- 
shift intervals. The values are plotted at the central redshift of the 
intervals. The dashed lines are the best fit densities of the LDDE 
model with an evolving JVh distribution depending on Lx an d Z 
(fit #4 in Tabled) . Data have been plotted using the " N ohs /7V mdl 
method" (see §4.1). 
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the fraction of absorbed objects with luminosity, and an 
increase with redshift. Both TVh distributions in which 
we allowed for a dependence of the absorbed objects on 
redshift or luminosity only (fits #2 and #3) are rejected 
at more than 99.5% confidence level. On the contrary, 
model #4 (Figures |S] and |5J) , where both a dependence 
on redshift and luminosity is allowed (see Figures ITU1 and 
, provides a very good representation of the data with 
a Xn h probability of 83%. 

4.1.1. Analysis of the uncertainties and systematic 

biases 

We analyzed how much our results could be affected 
by uncertainties in the completeness correction method 
used. These uncertainties could be introduced by errors 
in the spectroscopic classification of the AGN, and by 
the assumption that the measured fraction of AGN1 as a 
function of Lx and Nf/o and the derived two Lx~Lr re- 
lationships for AGN1 and AGN2 (see §3.3) hold also for 
the higher redshift, optically fainter, unidentified popu- 
lation. In order to measure the maximum allowed range 
of the HXLF parameters due to uncertainties in these 
assumptions, we have carried out the HXLF fits under 
the two very extreme hypotheses that all the unidentified 
AGN would follow either the Lx~Lr relationship typical 
of the AGN1 (eq. HJ, or the Lx-Lr, relationship typical 
of the AGN2 (eq. EJ). It resulted that the best-fit param- 
eters changed within the measured la uncertainties. The 
results are shown in Figure 03 where the largest allowed 
AGN density regions due to the uncertainties introduced 
by the completeness correction method used arc shown 
by dashed areas. 

About 60% of the AGN used in our analysis have their 



Fig. 10. — Nn distributions in various luminosity and redshift 
intervals. Top. High luminosity and low luminosity AGN at z<l. 
Bottom. High and low redshift AGN. The continuous lines are the 
observed distributions, the dotted lines are the assumed JVjj distri- 
butions (fit #4 in Table|2] evolving Nn distribution with a LDDE 
HXLF evolution), while the dashed lines are the expectations tak- 
ing into account the selection effects. 




42 43 44 45 46 0.5 1 1.5 2 2.5 3 3.5 4 
Log Lx (2-10 keV) z 

Fig. 11. — Observed fraction of absorbed (A r H >10 22 cm -2 ) 
AGN as a function of Lx and z. The dotted lines are examples 
of the intrinsic assumed distributions at various luminosities and 
redshifts (LDDE model, fit #4 in Table0. The long dashed lines 
are the corresponding average intrinsic assumed distributions of the 
sample used. The short dashed lines are the expectations taking 
into account the selection effects. 
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Ah column densities derived from the hardness ratios 
(those belonging to the AMSSn, H2XMM0.5, CDF-N 
and CDF-S samples). This method could introduce some 
systematic bias. For example, our simple absorbed power 
law model could tcndcntially underestimate the real col- 
umn densities, because scattered X-rays and circum- 
nuclear starburst X-rays can provide additional flux. 
This effect is expected to be stronger at lower luminosi- 
ties were the fraction of the light coming directly from 
the nucleus should be smaller. If this is the case, the ob- 
served decrease of the fraction of absorbed AGN with the 
intrinsic luminosity should be even stronger. It should 
be noted, however, that Perola et al. (2004) found a 
rather satisfactory correlation between the column den- 
sities measured from the hardness ratios and from the 
X-ray spectral fits in the HELLAS2XMM sample. 

Recently Tozzi et al. (2005) have published Ah mea- 
surements on the CDF-S sample, obtained using X- 
ray spectral fits. We took advantage of these measures 
to check whether the hardness ratio method introduces 
some relevant systematic bias. No relevant difference 
or systematic trend on either luminosity or redshift was 
found. In a subsample of z <1.2 AGN, using the hard- 
ness ratios we measure a fraction of 17/32 absorbed AGN 

143 



with L x >10 43 erg s"\ while Tozzi et al. (2005) find 
16/32. At lower luminosities (10 41 <Lx<10 43 erg s _1 ) we 
measure a fraction of 20/31 absorbed AGN (Lx >10 43 ), 
while Tozzi et al. (2005) find 18/31. In a subsample with 
10 43 < J L X <10 45 erg s _1 at redshift below 1.5 we measure 
a fraction 23/43 absorbed AGN, while Tozzi et al. (2005) 
find 24/43. At redshift above 1.5 we find a fraction 31/43 
absorbed AGN, while Tozzi et al. (2005) find 34/43. 

We also checked whether our results might depend on 
the assumed X-ray K-correction (see §3.2). We repeated 
the fit # 4 assuming T — 1.7 or T — 1.9, or assuming 
an exponential cutoff at energy Ec = 300 keV. It turned 
out that the changes of the parameters are within the la 
uncertainties. 

4.2. PLE model 

We also checked if a simpler pure luminosity evolution 
model were consistent with the data. By introducing the 
evolution factor 



e(z) 



(l + z)P 1 (z<z c ) 
e(z c ) (z > z c ), 

the PLE model is expressed as 



d$(L x ,z) d$(L x /e(z),0) 



dLogL> 



dLogLj 



(12) 



(13) 



The PLE fit (FigureU3and fit # 6 in Table© provides 
a less probable solution for the HXLF. Furthermore the 
PLE fit finds that the evolution stops at z c =1.08to;o6- 
This low value should be attributed to the fact that there 
is an increase with luminosity of the redshift peak of 
the density of AGN. Low luminosity (Lx<10 43 erg s~ 4 ) 
AGNs peak at z = 0.5, while high luminosity (Lx>10 46 
erg s _1 ) AGNs peak at z ~ 2. In this framework the 
PLE fit finds a weighted mean of the different redshift 
cut off values of the low and high luminosity AGNs. 

Although formally acceptable, z c =1.08io 06 ^ s signifi- 
cantly smaller than the previous estimates for the evo- 
lution of AGN1 in the hard X-rays (z c =2.4 ± 0.5; La 
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Fig. 12. — Density of AGN in luminosity bins as a function 
of redshift. The dashed lines are the best fit values of the PLE 
model with an evolving TVji distribution (fit #7 in Table 121 . Data 
have been plotted using the "7V obs /A rmdl method" (see §4.1). The 
dashed areas are the largest allowed regions due to uncertainties in 
the completeness correction method used (see §4.1.1). 



Franca et al. 2002), and in the optical band (z c ~2.0; see 
e.g. Boyle et al. 2000). This difference should be at- 
tributed to the fact that both the hard X-ray AGN1 and 
the optical QSO populate preferentially the bright part 
of the HXLF (see e.g. Figure and related discussion) 
which, also in the LDDE model, faces a redshift cut off 
larger than 1.5-2. If the fit of the PLE model is carried 
out with a fixed z c =2.0, it turns out unacceptable, with 
a xIf probability of 3.4xl0~ 6 % and a 2D-KS probabil- 
ity of 0.27%. On the basis of these results and on the 
fact that the PLE model over-predicts the 2-10 keV X- 
ray background (XRB, see Table 2) and the soft X-ray 
counts (see §4.3), we consider such a parameterization of 
the HXLF evolution to be ruled out. 

4.3. The Counts 

Down to the flux limit adopted in this analysis (10~ 15 
erg cm~ 2 s _1 ), the 2-10 keV counts predicted by the 
models described in Table El are in good agreement with 
both the counts of the whole sample (including the ob- 
jects without spectroscopic identifications), and with the 
Bauer et al. (2004) compilation (see Figure IT3*|l . The fit 
of Moretti et al. (2003) is also shown. This is an a pos- 
teriori test implying that our method, used to correct 
for the spectroscopic incompleteness of the faint sam- 
ples, is reliable. At faint fluxes (S'<10 _15 erg cm~ 2 s _1 , 
where there are no data in our samples) the LDDE model 
is consistent within the errors with the data, while the 
PLE model tends to over-predict the measured number 
density 12 . This is mainly due to a higher density of low 

12 As explained in §3.4, we have included in the Nn distribution 
a fraction of objects with 10 25 <A r H<10 26 cm -2 , equal to that in 
the interval 10 2 *<7V H <10 25 . 
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Fig. 13.— Differential counts of AGN in the 2-10 keV band. 
The filled circles are our estimates from the samples used in this 
analysis. The open squares are the estimates from Bauer ct al. 
(2004). The continuous line are the counts predicted by the LDDE 
model (fit # 4), while the short dashed line are the counts of the 
PLE model (fit # 6). The long dashed line are the estimates from 
a compilation of Moretti et al. (2003). 
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Fig. 15. — Predicted counts (from our best-fit LDDE model 
#4) of AGN for un-absorbed and absorbed AGN divided into two 
luminosity classes. The filled dot is the measure of the density 
of QS02 by Perola et al. (2004), the open circles are the density 
of QS02 derived by Padovani et al. (2004), while the triangle is 
the measure of the density of QS02 from the HBS28 sample of 
Caccianiga et al. (2004). 
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Fig. 14. — Differential counts of X-ray sources in the 0.5—2 
keV band (open squares) and of AGN-1 (filled squares) from a 
compilation of Hasinger et al. (2005). The continuous line are the 
counts predicted by the LDDE model (fit # 4), while the dashed 
line are the counts of the PLE model (fit # 6). 



luminosity AGNs in the HXLF (in comparison with the 
LDDE model) , and to the absence of a counter-evolution 
at high redshift in the PLE model. 
Although this analysis is based on observations made 



in the 2-10 keV band, it is instructive to compare our re- 
sults with the 0.5-2 keV counts. Of course, we should be 
aware that, when predicting the 0.5-2 keV counts, our re- 
sults depend on the spectral assumptions (a r=1.8 spec- 
tral slope plus photoelectric absorption) , which could not 
be valid below 2 keV. The PLE model over-predicts the 
observed soft counts as compiled by Hasinger ct al. (2005) 
at faint (5'<10 _15 erg cm~ 2 s _1 ) fluxes (see Figure ITU . 
The situation is even worse, since at faint fluxes we ex- 
pect a relevant contribution from normal X-ray galax- 
ies to the counts (about 20 deg -2 and 400 deg -2 at 
>5'o.5-2=10~ 15 erg cm -2 s _1 and S . 5 -2=10~ 16 erg cm' 2 
s _1 respectively; Ranalli et al. 2003, Bauer et al. 2004). 
On the contrary, the LDDE model provides a more ac- 
ceptable solution. At bright fluxes the observed counts 
are above our predictions because the X-ray sources are 
dominated by stars and clusters of galaxies which are not 
included in our models. 

4.3.1. The density of absorbed AGN 

The 2-10 keV predicted counts obtained from the best 
fit model for the HXLF (# 4) are shown in Figure ^1 
after being splitted according to X-ray absorption and 
X-ray luminosity. Most of the luminous (Lx>10 44 erg 
s _1 ), absorbed (Ah>10 22 cm -2 ) sources are AGN2 (see 
also the discussion in §3.3). Luminous, obscured AGN 
are usually referred to as QS02 and in the simplest ver- 
sion of the AGN unified scheme are predicted to be more 
numerous than QSOl by a factor comparable to that 
observed for lower luminosity Seyfert galaxies (about 3- 
4). Despite extensive searches, narrow-line optically lu- 
minous QS02 appear to be extremely rare and by far 
less numerous than broad line quasars (see Halpern, Er- 
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acleous & Forster 1998 and references therein). Because 
of the selection effects due to obscuration, X-ray surveys 
are expected to provide an unbiased census of the QS02 
population and a more reliable estimate of their space 
density. Several QS02 candidates (i.e. luminous, X-ray 
obscured sources) have been discovered by Chandra and 
XMM-Newton surveys (~ 20 in the HELLAS2XMM sur- 
vey, Fiore et al. 2003, Mignoli et al. 2004, Cocchia et 
al. 2005; - 30 in the CDF-S+CDF-N; a dozen in the 
CLASXS survey, Barger et al. 2005). A sizable fraction 
of them (from 50% to 75%) has been confirmed by deep 
optical spectroscopic observations. A notable example 
has been reported by Norman et al. (2002). It should be 
noted that the quality of spectroscopic observations is 
not uniform and, given the relatively high redshifts, the 
Ha and H/3 wavelengths are poorly sampled, thus ham- 
pering a "pure" optical classification. On the other hand, 
it is important to remind that the QS02 classification is 
wavelength dependent. Several, X-ray obscured, lumi- 
nous QS02 do not show any evidence of strong emission 
lines even in high quality optical spectra, among them 
the QS02 prototype NGC 6240 (Vignati et al. 1999). 
Keeping in mind these caveats, and adopting an admit- 
tedly arbitrary luminosity threshold (Lx>10 44 erg s _1 ), 
we obtain a QS02 space density of 60 and 267 deg -2 
at 2-10 keV fluxes brighter than 10 -14 and 10~ 15 erg 
cm~ 2 s _1 respectively. The slope of the integral counts 
significantly flattens below 10 -14 erg cm~ 2 s _1 . As a 
consequence, the QS02 surface density at fluxes much 
fainter than those actually probed (i.e. 10 -16 erg cm~ 2 
s _1 ) increases by a relatively small amount reaching 354 
deg -2 . These figures have at least a 5% error, corre- 
sponding to the HXLF normalization uncertainties. The 
predicted counts (Figure EJ are in good agreement with 
the QS02 space densities measured by Caccianiga et al. 
(2004), Perola et al. (2004) and Padovani et al. (2004). 

4.4. The XRB spectrum: a self- consistency check 

A detailed modeling of the XRB spectrum over the 
full ~2-400 keV range is beyond the scope of this paper. 
Nevertheless, it is important to check that the evolving 
X-ray luminosity function and the Nh distribution de- 
rived in the previous sections match the XRB flux at 
least in the 2-10 keV energy range. To this end it is use- 
ful to remind that the XRB intensity below ^10 keV, as 
measured by several imaging X-ray telescopes, is likely 
to be affected by systematic errors. In Figure El a com- 
pilation of XRB measurements is reported. The maxi- 
mum difference is of the order of 30% between the EPIC- 
pn flux reported by De Luca & Molendi (2004) and the 
HEA01-A2 measure of Marshall et al. (1980). According 
to a recent reanalysis of the HEA01-A2 data (Revnivt- 
sev et al. 2004) the 3-60 keV spectrum should be renor- 
malized upward by about 15%. The resulting 2-10 keV 
flux is 1.96±0.10 erg cm -2 s _1 deg -2 . The solid curve in 
Figure El represents the integrated AGN spectrum ob- 
tained with our best fit LDDE model for the HXLF (fit # 
4) with the redshift and luminosity dependent Ah func- 
tion 13 . Our predicted 2-10 keV flux of 1.81xl0~ n erg 
cm~ 2 s _1 deg" 2 corresponds to ^92% of the Revnivtsev 

13 As explained in §3.4, we have included in the Nh distribution 
a fraction of objects with 10 25 <7Vh<10 2b cm -2 , equal to that in 
the interval 10 2 *<7V H <10 25 . 
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Fig. 16. — Integrated AGN spectrum computed from the best fit 
LDDE model for the HXLF (fit # 4) with the redshift and luminos- 
ity dependent Nh function. The model predictions are compared 
with a selection of XRB spectral measurements over the broad 0.5- 
400 keV energy range. The regions enclosed within the bow— ties 
correspond to the XRB spectrum and associated errors as mea- 
sured by ROSAT (Georgantopoulos et al. 1996) in the 0.5-2 keV 
band and by XMM-Newton (De Luca k, Molendi 2004) in the 2- 
8 keV band. The data points in the 3—15 keV energy range are 
from RXTE (Revnitsev et al. 2003), while those in the 3-60 keV 
are from HEA01-A2 (Marshall et al. 1980). The bow-tie at high 
energies (100-400 keV) is from HEAOl-Ki (Kinzer et al. 1997). 



et al. (2004) value and 108% of the original HEAOl- 
A2 measure. Given that the XRB synthesis has been 
obtained with very simple prescriptions for the intrinsic 
(before absorption) spectral energy distribution (a power 
law spectrum with r=1.8 plus an exponential high en- 
ergy cut-off e - E / Ec with E c = 200 keV for all AGN) it is 
reassuring to obtain a reasonably good description of the 
XRB spectral intensity over a broad energy range. As a 
final remark we note that an increasing ratio between ab- 
sorbed and un-absorbed AGN towards high redshifts has 
been already included in the synthesis models of Pompilio 
et al. (2000) and Gilli et al. (2001) though with different 
prescription for the absorption distribution. 

4.5. The LF of absorbed and un-absorbed AGN 

It is interesting to plot the evolution of absorbed 
(-/Vh>10 22 cm -2 ) and un-absorbed AGN, according to 
our best fit LDDE solution (fit # 4; Figure [T7J). As ex- 
pected, the absorbed AGN outnumber the un-absorbed 
ones at low luminosities and high redshifts. In the bot- 
tom panel of Figure El the HXLF is compared with the 
estimate of Miyaji et al. (2000) of the soft X-ray (0.5- 
2 keV) AGN LF (a slope T=1.8 has been assumed to 
convert the 0.5-2 keV luminosities into the 2-10 keV 
band). It turns out that, at low redshifts (z~0.25), the 
soft X-ray LF is almost coincident with our measure of 
un-absorbed AGN HXLF, while at high redshifts the soft 
X-ray LF is consistent, within the uncertainties, with the 
total HXLF. This behavior is explained by the stronger 
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Fig. 17. — Top. Density of absorbed and un-absorbed AGN in 
luminosity bins as a function of redshift (fit # 4). Bottom. Density 
of absorbed and un-absorbed AGN as a function of luminosity in 
four redshift intervals (fit # 4). The dot-dashed line is the 0.5-2 
keV LF of AGN of Miyaji et al. (2000; plotted assuming r=1.8). 

effects of absorption in the soft X-rays in comparison 
to hard X-rays, especially at low redshifts. As a conse- 
quence, at low redshifts, only un-absorbed AGN are de- 
tected in the soft X-ray band, while at high redshifts the 
bias reduces, and the same population which is observed 
in the 2-10 keV band is detected. As a consequence, the 
soft X-ray LF faces a stronger (LDDE) evolution than 
observed for the HXLF (see Hasinger et al. 2005). 

4.6. The LF of AGN 1 and AGN2 

The space density and evolution of AGN1 and AGN2 
can be estimated using the above described method. In 
order to correct for the spectroscopic incompleteness of 
the faint samples we had to compute the completeness 
function g(Lx, z, Nu, R) (see §3.3) which is based on the 
estimate of the probability of an AGN to appear as an 
AGN1 as a function of L%, Ah and z: Ql(Lx, z, Njj) 
(shown in Figure EJ. With this estimate in hand (and 
keeping in mind the uncertainties on the AGN1-AGN2 
optical classification discussed in §3.3), we can derive the 
AGN1 luminosity function: 

*i(ix,«)= / *(Lx,z)f(L x ,z;N B )Ql(L x ,N s ,z)dN s , 

(14) 

where, as discussed in §3.4, f(Lx,z;Nn) is the Ah dis- 
tribution. The AGN2 density can be derived by sub- 



FlG. 18. — Top. Evolution of the luminosity function of AGN1 
up to z = 3.5. The data are the estimates of the density of AGN1 
from La Franca et al. (2002). Bottom. Evolution of the luminosity 
function of AGN2 up to z = 3.5. The lines have the same meaning 
as in the top panel. The small picture shows the difference among 
the luminosity functions of AGN1 and AGN2 at z = 0.1 and z = 
2.5. 



stituting, in the above formula (eq. I14fl . Ql(Lx, Ah, z) 
with l-Ql(Lx, Ah, z). As can be seen in Figure ITB1 at 
low redshifts (z<0.5), the AGN2 density low luminosities 
(Lx r °10 42 -10 43 erg s _1 ) is about five times larger than 
that of AGN1, while the latter outnumber the former by 
an order of magnitude at high luminosities (Lx^lO 46 ). 
In Figure [H the AGN1 density in the 2-10 keV band 
from the BeppoSAX HELLAS survey as computed by La 
Franca et al. (2002) is reported. The present estimate 
of the AGN1 luminosity function is consistent, within 
the uncertainties, with the La Franca et al. (2002) find- 
ings, thus confirming that Ql is a reliable measure of the 
probability of an AGN to appear as an AGN1. 

5. DISCUSSION 

5.1. Comparison with previous results 

A specific procedure to correct for the spectroscopic 
incompleteness of faint X-ray sources which also takes 
into account the selection effects due to X-ray absorption 
has allowed us to use a large AGN sample to compute the 
HXLF. Our results extend those of Cowie et al. (2003) 
and Barger et al. (2005), where no correction for X-ray 
absorption is adopted, and the upper limits to the AGN 
density arc estimated by assigning to the unidentified 
sources the redshifts corresponding to the centers of each 



14 



LA FRANCA ET AL. 




10-" io- 13 io-' 2 10- 

S (2-10 keV) [erg s"' cm"'] 

Fig. 19. — Observed fraction of absorbed AGNs as measured 
from the samples used in this analysis. The lines correspond to the 
predictions of the fitted HXLF listed in Tablel2l The open circle is 
the value obtained from the sample of Grossan (1992; not included 
in the analysis of the HXLF). 



Lx~z bin. 

A LDDE model provides the best fit to the HXLF evo- 
lution up to z = 4, in agreement with the Ueda et al. 
(2003) findings obtained using a smaller and brighter 
sample, and also with the estimates of Cowie et al. 
(2003), Fiore et al. (2003), Hasinger et al. (2005) and 
Silverman et al. (2005) who found that the AGN number 
density for luminosities lower than ~ 10 44 erg s~ 1 peaks 
at lower redshifts than that of higher luminosity AGN. 

The new result or our analysis concerns the luminos- 
ity and redshift dependence of the fraction of absorbed 
AGN, which decreases with luminosity and increases 
with redshift. The luminosity trend was already pointed 
out by Ueda et al. (2003) (see also Hasinger et al. 2005). 
In the Ueda et al. (2003) best fit model the fraction of 
absorbed AGN is 57% and 36% at i x =10 42 5 erg s^and 



L x =10 45 erg s _1 respectively. Taking into account only a 
luminosity dependence on the fraction of absorbed AGN 
in our sample (fit # 2 in Table 2), the corresponding 
fractions at L x =10 42 5 erg s _1 and L x =10 45 erg s _1 are 
68% and 40% respectively. The two results are remark- 
ably similar, especially if we note that absorbed AGN 
with 24<LogAn<25 are included in our sample but not 
in Ueda et al. (2003). 

The increase of the fraction of absorbed AGN with 
redshift, instead, emerges only with our analysis. The 
difference with respect to the Ueda et al. (2003) findings 
is due to the larger sample extending to fainter fluxes 
used in the present analysis. Indeed, if we restrict our 
analysis to a subsample (fit #5 in Table 2) of 207 objects 
from the Piccinotti, AMSSn and CDF-N catalogs and 
thus quantitatively similar to that used by Ueda et al. 
(2003) the uncertainties become so large that the redshift 
dependence is no longer significant while the luminosity 
dependence is recovered. 



It is worth noting that the luminosity and redshift de- 
pendence of the absorbed AGN fraction would disappear 
if one flux limited sample only were analyzed (as dis- 
cussed by Perola et al. 2004). A flux limited sample 
selects low luminosity AGN at low redshifts (which, ac- 
cording to our analysis, are more absorbed) and high 
luminosity AGN at high redshift (which are more ab- 
sorbed as well!). Then the average fraction of absorbed 
AGN turns out to be roughly constant. Only combin- 
ing several samples, and thus covering wide strips of the 
Lx~z plane with almost constant redshift or luminosity, 
it is possible to disentangle the true dependencies. 

A simple AGN model based on the unified paradigm 
has been adopted by Treister et al. (2004). Assuming 
that obscured AGN outnumber unobscured ones by a 
factor 3 without any luminosity and/or redshift depen- 
dence they claim to be able to reproduce the observed 
counts, redshift and Ah distributions in the CDF-N and 
CDF-S samples once all the selection effects are properly 
taken into account. More recently Treister & Urry (2005) 
revised their previous analysis including a luminosity de- 
pendence of the fraction of absorbed AGN which appears 
to provide an equally good fit to several observational 
constraints. However, in both works, no comparison be- 
tween the predicted and observed Ah distributions as a 
function of both L x and z is made. We have repeated 
our analysis assuming the Treister & Urry (2005) Ah 
distribution (see their Fig. 1). Either using the CDF-N 
plus CDF-S samples only, or the full AGN sample used in 
this work, we found that the only statistically acceptable 
models are those including a dependence of the fraction 
of absorbed AGN as a function of I/ X and z, with a be- 
haviour similar to what measured in the present paper 
(see §4). 

The observed and predicted fractions of absorbed 
(A H >10 22 /Total) AGN as a function of the observed flux 
are shown in Figure 1191 The open circle is from the 
Grossan sample (1992; not included in this analysis 14 ), 
and it is plotted in order to show the uncertainties at 
bright fluxes. As already described by several authors 
(e.g. Comastri et al. 2001; Tozzi et al. 2005; Perola et 
al. 2004 and references therein) the average X-ray spec- 
trum significantly hardens towards faint fluxes and this 
change is mostly concentrated in the 10~ 14 -10~ 13 erg 
cm~ 2 s _1 range, where the fraction of absorbed AGNs 
rises from about 20% to about 50%. For this reason, this 
measure is a very powerful tool to discriminate between 
different evolutionary scenarios for the Ah distribution. 
The only acceptable description of the observed ratio be- 
tween absorbed and un-absorbed AGN as a function of 
the hard X-ray flux is obtained only if the ratio depends 
on both luminosity and redshift (fits #4 and #6 in table 
2 for LDDE and PLE, respectively). These two models 
are indistinguishable, and in fact the fraction of absorbed 
AGN as a function of the X-ray flux is less sensitive to 
the shape and evolution of the HXLF than to the evolu- 
tion of the Ah distribution (see also §3.5). 

Recently, Alexander et al. (2005a, 2005b) have found 
evidences that a fraction of the z>l submillimiter cmit- 

14 A proper reassesment of the Grossan sample seems necessary 
before using it extensively for a detailed statistical analysis. Sec 
e.g. Bianchi et al. (2005) for a discussion on a few sources of the 
sample. 
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Fig. 20. — Top. Intrinsic (before absorption) luminosity density 
in the 2-10 keV band as a function of redshift from our best fit 
HXLF. Bottom. Total accreted mass as a function of redshift. 



ting galaxies harbor obscured AGN. They argued that 
the black holes are almost continuously growing through- 
out vigorous star-formation episodes. These results 
are in agreement with the hydrodinamical simulation of 
galaxy mergers by Di Matteo et al. (2005) and Springel 
et al. (2005), where the growth of both the black holes 
and stellar components are taken into account. In this 
framework, our result of an increase of the fraction of 
absorbed AGN with the redshift is in agreement with a 
picture where the peak epoch of the star formation (z=l- 
2) corresponds to a heavily obscured rapid black-hole 
phase, which is ultimately proceeded by an unobscured 
quasar phase (Alexander et al. 2005a, 2005b, Hopkins et 
al. 2005). 

5.2. Accretion history of the Universe 

Our measure of the HXLF cosmological evolution di- 
rectly constrains the history of the formation of super- 
massive black holes (SMBH) in the galactic centers not 
only for luminous un-obscured AGN, that can be traced 
also at longer wavelengths (optical, soft X-rays) , but also 
for the less luminous or obscured AGN. Starting from our 
best fit HXLF it is possible to derive the intrinsic (i.e. 
before absorption) luminosity density in the 2-10 keV 
band in the Universe as a function of redshift: 



J i x $(ix,z)dLogi X - 



(15) 



This quantity can be converted into the energy den- 
sity production rate per comoving volume by means of a 
bolometric correction factor K (Lbd = KL X ). The mass 
inflow rate onto a SMBH, M # , is related to the bolometric 
luminosity of the AGN, L^oi, by M. = Lboi{l — e )/ ec - 2 , 
where e is the radiative efficiency of the accretion flow 
(typically taken to be about 0.1; see e.g. Yu & Tremaine 



2002, Marconi et al. 2004, Barger et al. 2005). Once 
a value for e and K is adopted it is straightforward to 
derive the accretion rate density as a function of redshift: 



p.(z) 



1 



KL x 1>(L x ,z)dLogL x , (16) 



and the total accreted mass, i.e. the total density in 
massive black holes, if we assume that the initial mass of 
seeds black holes at redshift z s is negligible with respect 
to the total mass: 



Pbh(z) 



p.(z)-^dz. 

az 



(17) 



In the top panel of Figure |5D] we show our direct esti- 
mate of the intrinsic luminosity density in the 2-10 keV 
band as a function of redshift from our best fit HXLF. 
Assuming e = 0.1 and the luminosity dependent bolo- 
metric correction extensively discussed by Marconi et 
al. (2004; their eq. 21), the total density of massive 
black holes as a function of redshift reported in Fig- 
ure |201 is obtained by integration of the HXLF start- 
ing from z s = 4.5, for L x >10 41 erg s _1 and A H <10 26 
cm" 2 (as explained in §3.4). The final accreted mass 
correspond to a black hole mass density in the local Uni- 
verse of pbh = 3.2 hj x 10 5 M Mpc -3 . A somewhat 
higher value pbh — 4.0 h^ x 10 5 M Mpc" 3 is ob- 
tained for a single valued bolometric correction factor 
K = 40 (Elvis et al. 1994). These results are consistent, 
within the errors, with the SMBH density estimate of 
p BH — 4.6^2 4 /i 2 x 10 5 derived from dynamical studies 
of local galaxies bulges (see e.g. Marconi et al. 2004 and 
Ferrarese 2002). As shown in Figure QUI the vast ma- 
jority of the accretion rate density and black hole mass 
is produced by the low luminosity AGN (L x <10 44 -10 45 
erg s _1 ) down to redshift 2~1. As already shown by the 
LDDE model of the HXLF, high luminosity AGN are al- 
ready formed at redshift ~2 while low luminosity AGN 
keep forming down to 2~1. This result is in qualitative 
agreement with semi-analytical models for galaxy forma- 
tion and star formation rates, such as those of Balland et 
al. (2003), Menci et al. (2004), Granato et al. (2004), or 
with the hydrodinamical simulations such as those of Di 
Matteo et al. (2005), Springel et al. (2005), and Hopkins 
et al. (2005). 

6. CONCLUSIONS 

We have devised a method to compute the AGN HXLF 
which allow us to correct for both the spectroscopic in- 
completeness of the faint samples, and for the selection 
effects due to the X-ray K-correction. Thanks to this 
method we have been able to collect a sample of about 
500 AGN up to z —A. The most important results can 
be summarized as follows: 

• There is evidence that the fraction of absorbed 
(A r H >10 22 cm" 2 ) AGN decreases with the X-ray 
luminosity, and increases with the redshift. 

• the AGN HXLF up to z=4 is best represented by a 
LDDE model where the low luminosity (L x ~10 43 
erg s _1 ) AGN peak at z ~ 0.7 while high luminosity 
AGN (L x >10 45 erg s" 1 ) peak at z ~ 2. 
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• We can rule out a PLE model on the basis of several 
arguments which take into account the discrepan- 
cies with the optical and hard X-ray LF of AGN1, 
and the over-predictions of the soft X-ray counts 
and XRB intensity. 

• We estimate a density of supermassive black holes 
in the local Universe of pbh = 3.2 h^ x 10 5 Mq 
Mpc -3 , which is consistent with the recent esti- 
mates of local galaxies black hole mass function. 
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